Chaotic shock waves of a Bose-Einstein condensate 
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It is demonstrated that the well-known Smale-horseshoe chaos exists in the time evolution of 
the one-dimensional (ID) Bose-Einstein condensate (BEC) driven by the time-periodic harmonic 
or inverted-harmonic potential. A formally exact solution of the time-dependent Gross-Pitaevskii 
equation (GPE) is constructed, which describes the matter shock waves with chaotic or periodic 
amplitudes and phases. When the periodic driving is switched off and the number of condensed 
atoms is conserved, we obtained the exact stationary states and non-stationary states. The former 
contains the stable 'non-propagated' shock wave, and in the latter the shock wave alternately col- 
lapses and grows for the harmonic trapping or propagates with exponentially increased shock-front 
speed for the antitrapping. It is revealed that existence of chaos play a role for suppressing the blast 
of matter wave. The results suggest a method for preparing the exponentially accelerated BEC 
shock waves or the stable stationary states. 
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A shock wave in a compressible classical fluid is char- 
acterized by a steep jump in gas velocity and density due 
to the collisions of particles. The physics and applica- 
tions of the shock waves have been thoroughly investi- 
gated during the last century in different fields of physics 
HI) EL HI- Very recently, a different type of shock waves 
in a quantum fluid is explored 0, H, Hi S B 9 > where 
the quantum gas is a Bose-Einstein condensate (BEC) 
governed by the Gross-Pitaevskii equation (GPE). Ex- 
periments and numerical simulations [l^, fill . fl2| have 
depicted a BEC that exhibits the traveling fronts with 
steep gradients of shock waves. The corresponding ana- 
lytical studies d, [H, [3, [TH have also shown that a shock 
wave could develop in the attractive or repulsive BECs. 
However, as a nonlinear Schrodinger equation the driven 
GPE allows existence of the recognizably chaotic behav- 
ior [l6|, [It], Hl| • Therefore, various physical behaviors of 
the BECs may be affected by the chaos, that leads to the 
chaotic BEC solitons Hil HH , chaotic atomic popu- 
lations [H, [H, [Hj] , chaotic quantum tunneling [25|, Hi| , 
chaotic Bogoliubov excitations [13] , and the chaotic BEC 
collapse [28[ . All the above theoretical works are based on 
the analytical or numerical approximations to the GPE, 
because of the nonintegrability of the system with exter- 
nal potentials. The aim of this paper is to present an 
exactly analytical evidence of the chaotic and oscillating 
shock waves in ID attractive or repulsive BECs driven by 
the time-periodic harmonic or inverted-harmonic poten- 
tials, containing the most studied harmonically trapped 
BEC with zero driving. 

It is well-known that classical Smale-horseshoe chaos 
can exist in a parametrically driven Duffing system with 
the cubic nonlinear and harmonic or inverted-harmonic 
force [29|, H(| ■ The harmonic potential is a most widely 
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used trapping potential in the investigations of BECs 
[3TI [HJ • The expulsive parabolic potential has also been 
applied to study the time evolution of quantum tun- 
neling [§j| , deterministic chaos [H, [Hj] and accelerated 
bright solitons [H, H3] (e.g. see Fig. 3B of Ref. [H 
for the attractive interaction) in BECs. The harmonic 
and inverted-harmonic potentials can be time-dependent 
[H, Hi! S3] ■ The presences of such potentials will break 
the integrability of GPE and bring the difficulties to the 
exact researches of BECs. Although for different initial 
and boundary conditions the nonlinear GPE have many 
solutions, the exact solution of GPE with harmonic or 
inverted-harmonic potential has not been reported yet 
even for the non-driving cases. Particularly, we know 
that a turbulent flow is a fluid regime characterized by 
chaotic, stochastic property changes, and as an impor- 
tant physical behavior the unstable and irregular quan- 
tum turbulence has been studied [4l[. Therefore, the 
chaotic shock wave as a different type of turbulent flow 
warrants further investigation. 

In the present paper, we treat the time evolution of a 
quasi-ID BEC created initially in a range near the poten- 
tial center and driven by the time-periodic harmonic or 
inverted-harmonic potential, and seek the formally exact 
solution of the time-dependent GPE. By using the exact 
solution we reveal that the well-known Smale-horseshoe 
chaos exists in the matter shock waves with unpredictable 
amplitudes, phases and wave fronts. When we switch off 
the periodic driving and keep the number of condensed 
atoms, the exact solution becomes explicit functions of 
spatiotemporal coordinates which contain the stationary 
and non-stationary states. The stationary states include 
the exact stable 'non-propagated' shock wave. The non- 
stationary states describe the exactly controllable shock 
waves which alternately collapse and grow for the har- 
monic trapping and repulsive interaction or propagate 
with exponentially increased shock-front speed for the 
antitrapping and attractive interaction. The suppression 
of chaos to the blast of matter wave is revealed. The re- 
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suits can be observed experimentally [3.ITH. and supply a 
method for preparing the stable stationary state and the 
exponentially accelerated shock waves which are similar 
to the accelerated bright soliton [36| . 

There are some different methods for the reductions 
from original 3D GPE to quasi- ID GPE. The result- 
ing equations contain the nonpolynomial version [42| 
and the cubic nonlinear one k0| which has been used 
in the regime of shock wave [JM Q3- Assuming the 
transverse wave function to be the ground state of a 
harmonic oscillator of frequency u> r , the governing lon- 
gitudinal GPE becomes the cubic nonlinear equation 

ihip t = -7^ip xx + \Ya(x,t) + g' ld \^\ 2 ]ip, where m is 
atomic mass, g' ld — mio r g§l '(2irK) = 2hio r a s is the quasi- 
1D atom-atom interaction intensity with a s being the s- 
wave scattering length, V^(x,t) — {a\m,J 1 x J r V[ coswt)x 2 
denotes the harmonic (a = 1) and inverted- harmonic 
(a = —1) potentials of frequency lo x with the driving 
of strength V{ and frequency u> [40j ■ Taking an experi- 
mentally suitable frequency ujo as the units of frequencies 
uj x , ui and normalizing the time, space and wave function 
with Wq , Iq = \Jhj (muj ) and l/Vki, the GPE becomes 
the dimcnsionless one 

#t = ~\^xx + [k a {t)x 2 + gi d \if)\ 2 ]ip, 

k a {t) = (a^-u>l + V 1 cos ut). (1) 

Here the interaction is reduced to gi c i = 2uj r a s j '(ljqIq) 
and the potential strength V\ is normalized by frio®. 
Throughout the paper we take \g\ d \ = 0.4 which means 
Hujq = 25muj 2 a 2 with m, ui r and a s determined by ex- 
periment. 

Noticing that Eq. (1) can describe a symmetry BEC 
system in the transformation x — > —x, we consider a 
real physical process in which a BEC is created exper- 
imentally with the symmetry profile between positions 
±£0 — at initial time to, then the condensed atoms 
propagate along ±x directions such that the bound- 
ary coordinates ±xo(i) are time-dependent. Letting the 
number of condensed atoms conserve as N, we face the 
definite problem with the initial data [f| 

fl>(x,to) = (t iX,to) I " l'1cr' (2) 
YK ' uy 10 for \x\ > L w 

and the boundary-dependent normalization condition 
1-1%) \ip(x,t)\ 2 dx = N for any t. For a practical BEC 
the boundary density \if}(±L, t o)\ 2 could be made nonzero 
by using Feshbach resonances [T4I that leads to the steep 
jumps in atomic density and flow velocity, which charac- 
terize the feature of shock wave @ with the shock front 
±xo(t) at which the front gradients ip x (±xo , t) is dis- 
continuous. Although the initial data of ip(±L,to) and 
ip x (±L,to) cannot be measured accurately, in mathemat- 
ics, they are corresponded with an unique solution of 
Eq. (1) strictly. Hence, if we construct an exact solution 
ip(x,t) which obeys Eq. (2) and the normalization condi- 
tion experimentally, it will describe the exact shock wave 



of the system. For other initial and boundary conditions, 
of course, one could obtain different type of solutions. In 
fact, the accelerated bright soliton of the system has been 
reported for the initial soliton state [36| . 

We now give the formally exact solution of Eq. (1) as 

1>{x,t) - [A(t)x + iB(t)]e l ^ +b ^ x2 \ (3) 

where A(t), B(t), a(t), b(t) are real functions of time 
and obey the coupled equations 

a + gidB 2 = 0, b + 2b 2 = -g ld A 2 ~ k a (t) , 

A + 3Ab = 0, B + Bb = 0. (4) 

The solution (3) can be easily proved by inserting Eqs. 
(3) and (4) into Eq. (1) directly. The atomic density and 
flow velocity are associated with the macroscopic wave 
function ijj(x,t) = J~p e 10 through 

p(x,t) = |V>(x,i)| 2 = A 2 x 2 + B 2 , 

v(x,t) = *e x = -(2bx-—). (5) 

In order to evidence the solution (3), we have to solve 
Eq. (4) as follows. From the third equation of Eq. (4) 
we get b = -A/(3A), b = A 2 /(3A 2 ) - A/(3A). Applying 
these formulas to the second one of Eq. (4) produces the 
decoupled equation 

A = 3k a (t)A + 3g ld A 3 + — . (6) 

This is just the parametrically driven Duffing equation 
with a 'quadratic damping' term. It is well known that 
according to the Melnikov chaos criterion for a certain pa- 
rameter region the Smale-horseshoe chaos exists in such 
a system [2^, |3(| • To confirm the existence of chaos nu- 
merically, we adopt the MATHEMATICA code 

T = 2ir /to; e[{Aneu>-, vnew-}] := {A[T],v[T]}/. Flatten 
[NDSolve[{A'[t] ==v[t],v f [t] == 5/3 * v[t} 2 /A[t] 
+3(0.5au;l + V!Cos[Ljt})A[t} + 3g ld A[t} 3 , A[0] == Anew, 
v[0] == vnew}, {A, v}, {t, 0, T}}]; Do[pici = ListPlot 
[Drop[Nestlist[e, {Random[Real , {0.2, 1.5}], 
Random[Real, {-0.8, 0.8}]}, 6010], 10], {i, 1, 20}] 

with the parameters g\ d = —0.4, 3V\ — 20.2 and random 
initial conditions {A(0) G [0.2, 1.5], 1(0) G [-0.8, 0.8]} 
to plot two groups of Poincare sections on the 'phase 
space' (A, A) respectively for (a) aui 2 — 0.8, ui = 6.0 
and (b) aui 2 — —0.8, uj — 5.4. We find some periodic and 
chaotic orbits from each group consisting of 20 Poincare 
sections. Here any Poincare section is a discrete set of the 
phase space points at every period of the external poten- 
tial. In Fig. la and Fig. lb, we show the superpositions 
of 20 Poincare sections for each group respectively, where 
Fig. la with a = 1 means the harmonic potential case 
and Fig. lb with a = — 1 the case of inverted-harmonic 
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potential. In the both cases, if we change g-yj, < to 
9id > 0, the orbits will tend to infinities for some pa- 
rameter sets or tend to zero for some other parameter 
sets. The chaotic region of parameter space has not been 
found for the repulsive interaction case yet. Given so- 
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FIG. 1: The Poincare sections on the dimensionless 'phase 
space' (^4, ^4) for the parameters gid = —0.4, 3Vi = 20.2 
and (a) aw| = 0.8, w = 6.0, (b) aw* = -0.8, w = 5.4. 
Each figure consists of 20 Poincare sections associated with 
different initial conditions. 

lution A(t) from Eq. (6), the functions B(t), a(t), b(t) 
can be derived from Eq. (4) immediately. Inserting them 
into the exact solution of Eq. (3), the atomic density and 
phase gradient of Eq. (5) become explicit forms, which 
govern the spatiotemporal evolutions of the condensed 
atoms. The evolutions can be periodic for A(t) on the 
periodic orbits of Fig. 1 and can be chaotic for A(t) on 
the chaotic orbits of Fig. 1. 

For the considered BEC created initially near the po- 
tential center, the condensed atoms will propagate as a 
shock wave with time-dependent wave front coordinate 
xq (t) . At any hxed time the atomic density profile of Eq. 
(5) is a parabola with the steep gradients of shock fronts. 
The corresponding normalization condition reads 



N = 



° \^(x,t)\ 2 dx = \A 2 (t)xl + 2B 2 {t)x . (7) 

-x 



the relation A = AoB 3 with constant Aq such that Eq. 
(7) infers the formula 2A 2 [B 2 (t)x ] 3 + 6[B 2 {t)x ] = 3N, 
namely B 2 (t)xo — [A(t)/Ao] 2 ^ 3 xo = A is a constant and 
the coordinate Xo(t) reads 



x (t) = X[A /A(t)} 



2/3 



(8) 



where A is fixed by the algebraic equation 2A 2 X 3 + 6A = 
3N and Aq is determined by the initial condition xo(to) = 
L in Eq. (2). It is interesting noting that function A(t) 
in Eq. (8) may be periodic or chaotic so that the shock 
wave front Xo(t) may periodically or chaotically oscillate. 
For the periodic solution A(t), the driven BEC breathes 
in the spatiotemporal evolutions. In the chaotic solution 
case, the shock wave front Xo(t), BEC density p(x, t) and 
flow velocity v(x,t) become unpredictable. 

In the non-driving case with V± = and k a 
Eq. (6) becomes integrable. Applying the test for- 
mula A 2 = CiA 2 + C 2 A 4 to Eq. (6) and integrat- 
ing this equation, we obtain the undetermined constants 
C\ = — %k a , C2 = Qgid and the two exact solutions 
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A a (t) 



sech 



(9) 



respectively for g\d < 0, a = — 1 and gid > 0, a = 1, 
where to is an integration constant. The exact solutions 
of Eq. (9) can be directly proved by substituting them 
into Eq. (6) with constant k a . The A-i(t) of Eq. (9) 
denotes a sech-shaped solution for the attractive interac- 
tion and inverted-harmonic potential case. It is just the 
homoclinic solution which is associated with the Smale- 
horseshoe chaos of Melnikov criterion [2^, [3(| when the 
periodic perturbation appears. This agrees with the nu- 
merical result shown in Fig. 1. Noticing sech(zi) = seci, 
the A\ (t) means a sec-shaped one for the repulsive inter- 
action and harmonic potential case. Given Eq. (9), from 
Eq. (4) we easily write the undetermined functions 



f /A a \ 2 / 3 
a a (t) = -g xd J y-^-j dt, a = ±1 



(10) 



From the third and fourth equations of Eq. (4) we know 



with Aq being an integration constant adjusted by the 
initial conditions. Inserting Eq. (9) into Eq. (10), a sim- 
ple calculation can give the explicit forms of the functions 
B a (t), a a (t) and b a (t). Combining these functions with 
Eq. (3), we obtain the corresponding explicit forms of 
the exact non-stationary solutions immediately. 

Applying Eq. (9) to Eq. (8), we find that for the sech- 
shaped solution A-i(t) the shock wave front xo(t) is pro- 
portional to casb?/ 3 [xu) x (t— to)] which increases exponen- 
tially fast. However, for the sec-shaped solution A\ (t) the 
shock wave front Xo(t) is proportional to cos 2 ^ 3 [^uj x (t — 
to)] which oscillates periodically. The zero points of the 
periodic Xo(t) are just the singular points of A(t), where 
the BEC density p(x,t) = A\{t)x 2 + [A^t) / A w ] 2/3 be- 
comes infinity. By inserting Eq. (3) into Eq. (1) and 
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applying Eq. (4) to the resulting equation, we find that 
at the singular points of A(t) the original GPE (1) is 
fulfilled, since the infinite terms are offsetted each other. 
For a finite time the infinite density means the escape 
of solution (3) and describes the alternate collapses and 
growths of the blast matter wave physically [f| UJ. In the 
both cases, the traveling fronts have the steep gradients 
of shock waves [T5L [ill . , since normalization condition 
(7) implies that for a finite time ip[x = ±xo(£), t] — com- 
plex constants and ip[\x\ > \xo(t)\,t] — such that the 
front gradients tfj x [x = ±x (t),t] seem to be steep. No 
chaos exists for the non-driving case. 

Comparison between the driving and non-driving cases 
reveals that existence of chaos could play a role for sup- 
pressing the escape of solution and the blast of matter 
wave. In fact, from Fig. 1 we can see that in the driv- 
ing case with chaos, the function A(t) is finite and does 
not vanish for most of the orbits. These imply the corre- 
sponding shock front Xo(t) and density p(x, t) having no 
singular points. 

For the time-independent 'spring constant' k a , we have 
also found that Eq. (4) has the simple solutions 



a = -fit, 6 = 0, A 
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Id 



B 2 = (11) 
9U 



Substituting Eq. (11) into Eq. (3) yields the exact sta- 
tionary state wave functions 



ip a (x,t) 




(12) 



Here the harmonic potential with a = 1 is associated 
with the attractive interaction gid < and negative 
chemical potential \x, and the inverted-harmonic poten- 
tial with a = — 1 corresponds to the repulsive interac- 
tion gid > and positive chemical potential. The sec- 
ond spatial derivative of Eq. (12) vanishes that implies 
the zero kinetic energy of BEC and the 'non-propagated' 
shock wave with the steep front gradients. The con- 
stant front xq is determined by the algebraic equation 
—aui^Xg + 6/1x0 = SgidN from the normalization. Such a 
non-propagated shock wave can be realized for the time- 
independent potential V = a^uj 2 x 2 in the region \x\ < xo 
and V = u\uj 2 x x 2 ) outside [43|. Under such a finite po- 
tential, the stability of BEC can be determined by the 
known criterion. For a > and gid > 0, the stability 
criterion reads [44[ fi = fj, s = V(x ) + gid\tp(xo)\ 2 . This 
implies that for the harmonic potential and repulsive in- 



teraction the stationary state of Eq. (12) is stable, since 
the ip{xo) meets the stability criterion. 

In conclusion, for an atomic BEC created initially in 
a range near the potential center and driven by the 
time-periodic harmonic or inverted-harmonic potential, 
we have demonstrated that the classical Smale-horseshoe 
chaos certainly exists in the time evolutions of the sys- 
tem. The formally exact solution of the time-dependent 
GPE, whose amplitude and phase depend on the solu- 
tions of the famous Duffing equation with periodic driv- 
ing and quadratic damping, has been constructed, which 
describes the matter shock waves with chaotic or peri- 
odic amplitude and phase. When the periodic driving is 
switched off and the number of condensed atoms is con- 
served, we arrive at the most studied BEC systems with 
parabolic potentials. Then the exact solutions become 
the explicit functions of spatiotemporal coordinates and 
govern the exact non-stationary states or the exact sta- 
ble stationary states. The stationary states are called 
the 'non-propagated' shock wave. The non-stationary 
shock waves alternately collapse and grow for the har- 
monic trapping and repulsive interaction or propagate 
with exponentially increased shock-front speed for the 
antitrapping and attractive interaction. It is revealed 
that existence of chaos play a role for suppressing the 
blast of matter wave. The results can be observed ex- 
perimentally and suggest a useful method for preparing 
the exponentially accelerated matter shock wave or the 
stable stationary state of the BEC system. 

The well-known criterion for the onset of temporal 
chaos is the Melnikov criterion based on the nonlinear 
ordinary differential equations. In order to apply such a 
criterion to demonstrate the chaotic behaviors of Eq. (1), 
we have to consider the ansatz (3) governed by the re- 
duced ordinary differential equations (4) and (6). In gen- 
eral, the investigation of spatiotemporal chaos also can 
be performed directly from the nonlinear partial differen- 
tial equation (1), by using the Deissler-Kaneko criterion 
, which relies on the determination of the time evo- 
lution of a function defined by the integral of the square 
modulus of the difference between wave functions with 
nearby initial conditions [281 ]. 
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